library(Seurat)
library(fastTopics)
library(RcppParallel)
library(Matrix)
library(dplyr)
library(RColorBrewer)
library(ggplot2)
library(ggExtra)
library(cowplot)
library(wesanderson)
#Set ggplot theme as classic
theme_set(theme_classic())Hem.data <- readRDS("../QC.filtered.cells.RDS")DimPlot(Hem.data,
reduction = "spring",
cols = c(wes_palette("FantasticFox1"),"grey60"),
pt.size = 0.5) & NoAxes()# Extract apical progenitors
Progenitors.data <- subset(Hem.data, idents = c(0,1,3))
DimPlot(Progenitors.data,
reduction = "spring",
pt.size = 0.5,
cols = c(wes_palette("FantasticFox1")[c(1,2,4)]),
split.by = 'ident') + NoLegend() & NoAxes()rm(Hem.data) ; gc()## used (Mb) gc trigger (Mb) max used (Mb)
## Ncells 3236850 172.9 5989767 319.9 4708438 251.5
## Vcells 338565731 2583.1 1044947158 7972.4 1015766977 7749.7
For this analysis we will keep only genes detected in at least 20 over 12325 cells
progenitors.counts <- GetAssayData(object = Progenitors.data[["RNA"]], slot = "counts")
dim(progenitors.counts)## [1] 18268 12325
num.cells <- Matrix::rowSums(progenitors.counts > 0)
genes.use <- names(x = num.cells[which(x = num.cells >= 20)])
progenitors.counts <- progenitors.counts[genes.use, ]
dim(progenitors.counts)## [1] 15314 12325
gc()## used (Mb) gc trigger (Mb) max used (Mb)
## Ncells 3224155 172.2 5989767 319.9 4708438 251.5
## Vcells 414394458 3161.6 1044947158 7972.4 1015766977 7749.7
set.seed(1)
fit <- fit_topic_model(t(progenitors.counts),
k = 15,
numiter.main = 200,
numiter.refine = 200,
method.main = "em",
method.refine = "scd",
control.main = list(numiter = 4, nc= 6),
control.refine = list(numiter = 4, nc= 6, extrapolate = TRUE),
verbose = "progressbar")## Initializing factors using Topic SCORE algorithm.
## Initializing loadings by running 10 SCD updates.
## Fitting rank-15 Poisson NMF to 12325 x 15314 sparse matrix.
## Running 200 EM updates, without extrapolation (fastTopics 0.5-59).
## Refining model fit.
## Fitting rank-15 Poisson NMF to 12325 x 15314 sparse matrix.
## Running 200 SCD updates, with extrapolation (fastTopics 0.5-59).
# Add cells' topics loading to the metadata
Progenitors.data@meta.data <- cbind(Progenitors.data@meta.data, fit$L)FeaturePlot(object = Progenitors.data,
features = paste0("k", 1:15),
cols = rev(brewer.pal(10,"Spectral")),
reduction = "spring") & NoLegend() & NoAxes()FeaturePlot(object = Progenitors.data,
features = paste0("k", c(15,12,9,8,14,6)),
cols = rev(brewer.pal(10,"Spectral")),
reduction = "spring",
order = T) & NoLegend() & NoAxes()set.seed(1)
pca <- prcomp(fit$L[,c(15,12,9,8,14,6)])$x
clusters <- cluster::pam(pca, k = 6)$clusteringProgenitors.data@meta.data$TopicsKmeans <- as.numeric(clusters)
FeaturePlot(object = Progenitors.data,
features = "TopicsKmeans",
cols = c(wes_palette("FantasticFox1"),"grey90", "grey40"),
reduction = "spring") & NoLegend() & NoAxes()Idents(Progenitors.data) <- Progenitors.data$TopicsKmeans
DimPlot(Progenitors.data,
reduction = "spring",
pt.size = 0.5,
cols = c(wes_palette("FantasticFox1"),"grey90", "grey40"),
split.by = 'ident') + NoLegend() & NoAxes()domain.markers <- FindAllMarkers(Progenitors.data,
only.pos = F,
min.pct = 0.25,
logfc.threshold = 0.3)
domain.markers.filtered <- domain.markers %>%
group_by(cluster) %>%
top_n(n = 4, wt = avg_log2FC) %>%
unique()DotPlot(Progenitors.data,
features = unique(domain.markers.filtered$gene)
) + RotatedAxis()DotPlot(Progenitors.data,
features = c("Shisa2", "Wif1", "Rassf4", "Sulf1", "Dkk3", "Dlk1","Ttr")
) + RotatedAxis()ident = c("Dorso-Medial_pallium", "ChP", "Medial_pallium", "Hem", "ChP_progenitors", "Thalamic_eminence")
Progenitors.data$progenitor_type <- sapply(Progenitors.data$TopicsKmeans,
FUN = function(x) {x= ident[x]})
Idents(Progenitors.data) <- Progenitors.data$progenitor_typeDimPlot(Progenitors.data,
reduction = "spring",
pt.size = 0.5,
cols = c(wes_palette("FantasticFox1"),"grey90"),
split.by = 'ident') + NoLegend() & NoAxes()Hem.data <- readRDS("../QC.filtered.cells.RDS")Hem.data$Cell_ident <- sapply(Hem.data$Barcodes,
FUN = function(x) {
if (x %in% Progenitors.data$Barcodes) {
x = Progenitors.data@meta.data[x, "progenitor_type"]
} else {
x = paste0("seurat_clusters_", Hem.data@meta.data[x, "seurat_clusters"])
}
})DimPlot(object = Hem.data,
group.by = "Cell_ident",
reduction = "spring",
cols = c("#83c3b8", #"ChP"
"#009fda", #"ChP_progenitors"
"#68b041", #"Dorso-Medial_pallium"
"#e46b6b", #"Hem"
"#e3c148", #"Medial_pallium"
"#b7d174", #2
"grey40", #4
"black", #5
"#3e69ac" #"Thalamic_eminence"
))saveRDS(Hem.data, "../QC.filtered.cells.RDS")#date
format(Sys.time(), "%d %B, %Y, %H,%M")## [1] "19 novembre, 2021, 19,22"
#Packages used
sessionInfo()## R version 4.1.2 (2021-11-01)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.3 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0
##
## locale:
## [1] LC_CTYPE=fr_FR.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=fr_FR.UTF-8 LC_COLLATE=fr_FR.UTF-8
## [5] LC_MONETARY=fr_FR.UTF-8 LC_MESSAGES=fr_FR.UTF-8
## [7] LC_PAPER=fr_FR.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=fr_FR.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] wesanderson_0.3.6 cowplot_1.1.1 ggExtra_0.9 ggplot2_3.3.5
## [5] RColorBrewer_1.1-2 dplyr_1.0.7 Matrix_1.3-4 RcppParallel_5.1.4
## [9] fastTopics_0.5-59 SeuratObject_4.0.3 Seurat_4.0.5
##
## loaded via a namespace (and not attached):
## [1] plyr_1.8.6 igraph_1.2.8 lazyeval_0.2.2
## [4] splines_4.1.2 listenv_0.8.0 scattermore_0.7
## [7] digest_0.6.28 invgamma_1.1 foreach_1.5.1
## [10] htmltools_0.5.2 SQUAREM_2021.1 fansi_0.5.0
## [13] magrittr_2.0.1 tensor_1.5 cluster_2.1.2
## [16] ROCR_1.0-11 limma_3.50.0 recipes_0.1.17
## [19] globals_0.14.0 gower_0.2.2 matrixStats_0.61.0
## [22] MCMCpack_1.6-0 spatstat.sparse_2.0-0 prettyunits_1.1.1
## [25] colorspace_2.0-2 ggrepel_0.9.1 xfun_0.28
## [28] crayon_1.4.2 jsonlite_1.7.2 spatstat.data_2.1-0
## [31] survival_3.2-13 zoo_1.8-9 iterators_1.0.13
## [34] glue_1.5.0 polyclip_1.10-0 gtable_0.3.0
## [37] ipred_0.9-12 MatrixModels_0.5-0 leiden_0.3.9
## [40] future.apply_1.8.1 abind_1.4-5 SparseM_1.81
## [43] scales_1.1.1 miniUI_0.1.1.1 Rcpp_1.0.7
## [46] progress_1.2.2 viridisLite_0.4.0 xtable_1.8-4
## [49] reticulate_1.22 spatstat.core_2.3-1 stats4_4.1.2
## [52] lava_1.6.10 prodlim_2019.11.13 truncnorm_1.0-8
## [55] htmlwidgets_1.5.4 httr_1.4.2 ellipsis_0.3.2
## [58] ica_1.0-2 farver_2.1.0 pkgconfig_2.0.3
## [61] nnet_7.3-16 sass_0.4.0 uwot_0.1.10
## [64] deldir_1.0-6 utf8_1.2.2 caret_6.0-90
## [67] labeling_0.4.2 tidyselect_1.1.1 rlang_0.4.12
## [70] reshape2_1.4.4 later_1.3.0 munsell_0.5.0
## [73] tools_4.1.2 generics_0.1.1 ggridges_0.5.3
## [76] evaluate_0.14 stringr_1.4.0 fastmap_1.1.0
## [79] yaml_2.2.1 goftest_1.2-3 mcmc_0.9-7
## [82] ModelMetrics_1.2.2.2 knitr_1.36 fitdistrplus_1.1-6
## [85] purrr_0.3.4 RANN_2.6.1 pbapply_1.5-0
## [88] future_1.23.0 nlme_3.1-153 mime_0.12
## [91] quantreg_5.86 compiler_4.1.2 plotly_4.10.0
## [94] png_0.1-7 spatstat.utils_2.2-0 tibble_3.1.6
## [97] bslib_0.3.1 stringi_1.7.5 highr_0.9
## [100] lattice_0.20-45 vctrs_0.3.8 pillar_1.6.4
## [103] lifecycle_1.0.1 spatstat.geom_2.3-0 lmtest_0.9-39
## [106] jquerylib_0.1.4 RcppAnnoy_0.0.19 data.table_1.14.2
## [109] irlba_2.3.3 conquer_1.2.1 httpuv_1.6.3
## [112] patchwork_1.1.1 R6_2.5.1 promises_1.2.0.1
## [115] KernSmooth_2.23-20 gridExtra_2.3 parallelly_1.28.1
## [118] codetools_0.2-18 MASS_7.3-54 withr_2.4.2
## [121] sctransform_0.3.2 hms_1.1.1 mgcv_1.8-38
## [124] parallel_4.1.2 quadprog_1.5-8 grid_4.1.2
## [127] rpart_4.1-15 timeDate_3043.102 tidyr_1.1.4
## [130] coda_0.19-4 class_7.3-19 rmarkdown_2.11
## [133] ashr_2.2-47 Rtsne_0.15 mixsqp_0.3-43
## [136] pROC_1.18.0 shiny_1.7.1 lubridate_1.8.0
Institute of Psychiatry and Neuroscience of Paris, INSERM U1266, 75014, Paris, France, matthieu.moreau@inserm.fr↩︎